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PACS 64 . 70 . kp - Ionic crystals 
PACS 68 . 65 . Ac - Multilayers 

PACS 73 . 20 . -r - Electron states at surfaces and interfaces 

Abstract -We study the ground-state structures of identical classical point charges with Coulomb 
interactions, confined between two symmetric parallel charged walls. For the well understood ho- 
mogeneous dielectric case with no electrostatic images, the charges evenly condense on the opposite 
walls, thereby forming a bilayer Wigner crystal; five structures compete upon changing the inter- 
wall separation. Here, we consider a dielectric jump between the walls and a solvent in which 
charges are immersed, implying repulsive images. Using recently developed series representations 
of lattice sums for Coulomb law, we derive the complete phase diagram. In contrast to the ho- 
mogeneous dielectric case, the particles remain in a hexagonal Wigner monolayer up to a certain 
distance between the walls. Beyond this distance, a bifurcation occurs to a sequence of Wigner 
bilayers, each layer having a nonzero spacing from the nearest wall. Another new phenomenon is 
that the ground-state energy as a function of the wall separation exhibits a global minimum. 



> 

(N 

o 



X 



Introduction. — Bilayer Wigner crystals of electrons 
appear in various experimental settings: GaAs quantum 
wells [TJ[5] or other semiconductors [3J, quantum dots [3], 
dusty plasmas [5], colloids [5], etc. Theoretically, it is 
relevant and non-trivial to understand the ground-state 
features of Coulombic bilayers, starting with the classical 
limit. This subject in its homogeneous dielectric version 
has received significant attention in last years, as such 
[7^ll2j or within more general finite temperature analysis 

prog. 

We consider the ground state of an ensemble of (say 
elementary) classical point charges e, interacting through 
the 1/r usual Coulomb pair potential. The charges are 
confined between two symmetric parallel walls at distance 
d, carrying a uniform surface charge of density — ae (see 
fig. [TJ. The system as a whole is electroneutral. The 
static dielectric constant of the walls e' can be, in general, 
different from the constant e of the solvent in which the 
charges are immersed. The dielectric jump between the 
two moieties is defined as A = (e— e')/(e+e'), so that A G 
[— 1, 1]. The interval A € [—1, 0) corresponds to attractive 
and A G (0, 1] to repulsive electrostatic image charges. In 
realistic systems, in particular those of biological interest, 



the wall mimics the interior of a polarizablc colloid with 
e' < 10 while e ~ 80 for water solvent, so A is close to 1. 

In the well understood homogeneous dielectric case 
A = 0, Earnshaw theorem [T7] tells us that in the ground 
state, the charges evenly condense on the opposite walls. 
The structure of the formed bilayer depends on a single 
dimcnsionlcss parameter 77 = d^/a. For 77 = 0, a genuine 
hexagonal Wigner crystal (structure I) is stable (provides 
the true minimum of the energy) ; this structure becomes 
unstable for an arbitrary small n > |11I18019] . In the op- 
posite limit 77 — > 00, the two layers decouple and a hexag- 
onal Wigner crystal is formed on each plate (structure V); 
to minimize the inter-layer repulsion, these two crystals 
adopt a staggered configuration. For intermediate values 
of 77, three other structures are met, see fig. [5] a staggered 
rectangular lattice (structure II), a staggered square lat- 
tice (structure III) and a staggered rhombic lattice (struc- 
ture IV) [9]. The transitions between structures (II, III) 
and (III,IV) are of second order with mean-field critical 
indices [TMT§] . The transition between structures (IV,V) 
is discontinuous, characterized by a skip of the angle tp to 
7r/3 and by a mutual shift of the lattices on the opposite 
surfaces. 
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Fig. 1: The geometry of two dielectric walls at distance d, car- 
rying the surface charge density —ere. The particles of charge 
e form a Wigner bilayer which consists of two sub-lattices A 
and B at distance d; a monolayer corresponds to d = 0. The 
distance of each sub-lattice to the nearest wall is I = (d — d)/2. 
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Structures I, II and III 



Structure IV 



Structure V 



Fig. 2: Schematic representation of the different ground states 
encountered when the dimensionless distance r\ increases, in 
the homogeneous dielectric case A = 0. The open and filled 
symbols show the locations of ions on the opposite surfaces. 
The arrows are for the lattice vectors a\ and a?, from which 
we define the aspect ratio 8 — |a2|/|ai|: 5 = y3 for structure 
I, 5 = 1 for structure III, while structure II is intermediate, 
having 1 < 8 < ^/3. For structure IV, the order parameter 
is the angle <p between ai and 02. Structures IV and V have 
5 = 1. 



In this letter, we consider the dielectric jump Ae [0,1], 
implying repulsive electrostatic image charges. The pres- 
ence of repulsive image and self-image forces has intu- 
itively two important effects on the ground state: 

• the particles no longer collapse on the surface of the 
walls, but create a symmetric bilayer structure with 
a nonzero spacing from the walls, i.e. d ^ d in fig. [TJ 

• the intervals of r\, inside which the structures I-V are 
stable, are changed. In particular, a monolayer with 
structure I is expected to be stable up to a critical 
value 77 1— >n > beyond which it bifurcates into a 
bilayer with structure II. 

The nonzero spacing from the walls is confirmed in the 
limit 77 — > 00 where the two layers decouple (up to a stag- 
gered shift) and each of them forms a neutral entity with 
the corresponding charged wall. For the resulting one- 
wall geometry, it was shown in ref. |20j that, as a result 
of balance between attractive surface charge and repulsive 
image charge forces, the stable hexagonal Wigner crystal 
of particles is formed at some distance I from the wall. 
This distance is given by / = at/2, where a is the lattice 
spacing of the two-dimensional Wigner crystal determined 
by the electroneutrality condition as 1 = (v / 3/2)a 2 tr and 
the dimensionless parameter t is the solution of the tran- 
scendental equation 



{t 2 + j 2 + jk + fc 2 ) 3 / 2 v/3 A t' 



The lattice sum on the lhs can be expressed as 

- / du^le- ut2 [9 3 (e- u )9 3 (e~ 3u ) + f? 2 ( e -")0 2 (e- 3u )] , 
77 Jo 

(2) 

where 9 2 (q) = and 6 3 (q) = EZ-ool^ 

are the Jacobi theta functions with zero argument. 

Method. — Our goal is to calculate the Coulomb in- 
teraction energy of the electroneutral system pictured in 
fig. [TJ Each point will be represented by Cartesian coor- 
dinates as r = (R, z), where R = (x, y) and the z axis is 
perpendicular to the walls. The system is inhomogeneous 
along the z-axis, symmetric with respect to the z = 
plane. The wall surfaces, charged uniformly by — ae, are 
localized at z = ±d/2. The two Wigner sublayers, each 
formed of N/2 point charges e, are localized at z = d/2 
(sub-lattice A with sites {R^ 1 }) and at z = —d/2 (sub- 
lattice B with sites {R^}). The distance of the Wigner 
sub-lattice A (rcsp. B) to the nearest wall surface at 
z = +d/2 (resp. z = -d/2) is I = (d-d)/2. We have I > 
for the considered case with repulsive images A e (0, 1] 
and I = otherwise (for attractive images A € [—1,0), a 
small hard core around the particles is necessary to pre- 
vent thermodynamic collapse onto their images). 

The Coulomb energy consists of three contributions: 

(i) Particle-particle interactions: The direct Coulomb 
potential between two points r = (R, z) and r' = (R', z'), 
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situated in between the walls, is given by 
u (r,r ) 



e ^|R-R'| 2 + (z-z') 2 



(3) 



Due to the dielectric inhomogeneity, a unit charge at r = 
(R, z) has an infinite number of electrostatic images with 
the same vector R [3T|. The first sequence of images, 
generated by considering first the reflection with respect to 
the wall surface at z = +d/2, is localized at (— l) n+1 (nd — 
z) with associated charges A™ (n = 1, 2 . . .). The second 
sequence, generated by considering first the reflection with 
respect to the wall surface at z = —d/2, is localized at 
(— l) n (nd + z) with associated charges A™ (n = 1, 2 . . .). 
The image Coulomb potential thus reads 



A 



1 °° f 



-R'| 2 + (2nd + z - z') 2 
A 2 " 



v /|R-R'| 2 + (2n<i-z + z') 2 
A 2 ™ -1 

^/|R - R'| 2 + [(2n -l)d + z + z'} 2 

\2n-l 

v/|R - R'| 2 + [(2n -l)d-z- z'} 2 



(4) 

The notation Ui m (r, r') = Uj m (|R — R'|; z, z') will be often 
used in what follows. The total energy of particles and 
their images is given by 



^pp- 2 



j¥=k j,k 



(5) 



The term j = k in the second sum corresponds to the in- 
teraction of particle j with its own image. Let us choose 
as a reference particle the one at site 1 of lattice A lo- 
calized on the plane z = d/2. The energy per particle is 
expressible as 
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(6) 



Here, = |Rf -R/| and similarly i^/ 3 = |Rf -Rf |. 
Owing to the symmetry of the problem, the same formula 
holds if the reference particle sits on lattice B. 

The lattice sums in the above expression have to be 
"regularized" by background terms in the following sense. 



Let us consider two plates at distance z, each of a large 
surface S. They carry the fixed surface charge density 
— ere. At each plate, there are point charges e of surface 
density a. The total number of particles TV = 2a S, so 
that the system as a whole is electroncutral. The dielectric 
constant is homogeneous, equal to e. The energy per unit 
surface due to the background-background interaction of 
the two (—ere) plates at distance z is given by E^/S = 
—2Tr(ae) 2 z/e. Equivalently, E^/N = —nae 2 z/e. The 
interaction energy of charge e with (—ere) plate at distance 
z is 2irae 2 z/e. The total "background energy" per particle 
thus reads e(z) = Trae 2 z/e. Such term has to be added to 
and subtracted from each of the lattice sums in ((5]). The 
first term is unchanged, the second term becomes 



2e 
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{Rtff + d 2 



2ixad 



etc. The regularized sums (like the one in the square 
bracket) vanish in the large-distance limit. Interestingly, 
they are amenable to a special representation by quickly 
converging series expansions, with the aid of the recent 
method developed in 18,19 - Note that the rcgularization 
of image interactions u; m requires addition and subtrac- 
tion of an infinite number of background energies. After 
performing the rcgularization of all lattice sums in ([5]), we 
end up with 
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(7) 



where the star means "regularized" . 

(ii) Particle-surface charge interaction: We now 

study the direct and image interactions of charged par- 
ticles with the fixed surface charge densities on the walls. 
The surface charge —ae at z = —d/2 has one image 
A (—ere) at the same z = —d/2 and a series of images 
(1 + A)A™- 1 (-ere) at z = (-l)"(2n- l)d/2 (n = 2,3,...). 
The surface charge —ae at z = +d/2 has one image 
A(— ere) at the same z = +d/2 and a scries of images 
(1 + A)A n -\-ae) at z = (-l)"- 1 (2n - l)d/2 (n = 
2,3,...). After simple algebra, the interaction energy of 
one particle of charge e with both surface charges (—ae) 
and their images takes the form 



E ps 2-Kae 2 (1 + A) 2 



~N 



(1-A) 



(8) 



The same result is obtained, as it should be, when one 
computes the interaction energy of both surface charges 
(—ere) with one charge e and all its images, 
(iii) Surface charge-surface charge interaction: The 
energy per particle due to the direct and image interaction 
between surface charges turns out to be 



E ss 
N 



irae 2 (1 + A) 2 
e (1-A) 2 



d. 



(9) 
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Note that this value equals to —1/2 of E ps /N defined by 
©, in close analogy with the homogeneous A = case. 
For the total energy E = E pp + E ps + E ss , we get 

f = f + H£W m 

The infinite lattice sums in E* p /N ', which is the regular- 
ized version of © , can be expressed as quickly converging 
series of the special functions [TBI ITS] 

z v {x,y)= / -e-**e-^, y>0, (11) 

which represent a generalization to two-layer problems of 
the so-called Misra functions [22] used in single-layer lat- 
tice summations [23] • The cutoff of the series at the 5th 
term corresponds to the precision of more than 17 decimal 
digits; we keep this cutoff in the present work. As concerns 
the number of counted images in the image Coulomb po- 
tential (jl]), we cut the series at n — 5 as well; there exist 
various checks indicating that this cutoff preserves excel- 
lent accuracy of our results. 

We do not assume that the inclusion of repulsive images 
will lead to the appearance of a new structure, different 
from those presented in fig. [5] The energy per particle of 
a structure has to be minimized with respect to all struc- 
ture parameters. In comparison with the A = case, the 
bilayer distance d is an additional variational parameter. 
Among the five afore-mentioned structures, the one which 
provides the minimum of the energy E is chosen as stable. 

The appearance of three, four or more layered structures 
cannot be, in principle, excluded. Other approaches, like 
Monte Carlo (MC) simulations, are more suitable for de- 
tecting such complicated structures. We only know that 
the case A = and the two extreme limits d — > (struc- 
ture I) and d — > oo (structure V) for all A <G [0, 1] are 
treated correctly within our assumption: 

• As soon as A > 0, the stability along z-axis of a 
particle inside structure I is influenced "positively" 
by image/self-image effects and "negatively" by di- 
rect Coulomb interaction with all other particles. The 
image effects become stronger when the distance be- 
tween walls d decreases; there always exists a finite 
distance d when the image force overwhelms the di- 
rect repulsion by other particles. This is why struc- 
ture I is expected to be stable in a finite interval of d 
values, in contrast to the A = case. 

• In the opposite limit d — > oo, our energy per particle 
(fit)]) reduces to that of the hexagonal Wigner crystal 
with images for one-wall geometry. Indeed, in that 
limit, the following terms survive from the regularized 
version of ©: the sum l/Rff, A/[e{d - d)] 
from the n = 1 term of Ui m (0; d/2, d/2) and 

— V 1 + 2ira(d - d) 



from the n = 1 term of Yl&i u im(Rff\ d/2, d/2). 
Substituting these terms into Eq. (fTU)) and expressing 
everything in terms of the distance I of the hexagonal 
Wigner crystal from the wall, we recover the basic 
expression for the energy per particle (9) in |20j for 
the charge valence q = 1. The two Wigner crystals 
form, together with the corresponding charged wall, 
neutral entities shifted as in structure V to ensure the 
lowest energy. 

We shall work with dimensionless distances 

t? = dy/a, fj = dy/a, £ = l\fa, (12) 

such that £ = (r\ — fj)/2. In the asymptotic limit r\ — > oo, 
the dimensionless distance of the Wigner layer from the 
corresponding wall is given by 

e =»kn o-"* (13) 

where t is the A-dcpendcnt solution of Eqs. (UJ) and 
This formula provides a check of the accuracy of our results 
in the region of structure V as well as criterion for reaching 
the asymptotic one- wall regime. 




Fig. 3: The phase diagram in the (77, A) plane. The dashed 
line, lying entirely in the structure-I region, corresponds to the 
absolute ground state (minimum of the total energy). 



Numerical results. — We have performed numerical 
calculations using the Mathcmatica software. The phase 
diagram in the (77, A) plane is pictured in fig. [3] Gen- 
erally speaking, the transition values of rj exhibit a much 
stronger dependence on A than the ones of fj presented in 
table [JJ The relevant fact is that for every A > the par- 
ticles form the Wigner monolayer (fj = 0) with structure 
I up to some positive 771-^11 > 0, as was anticipated. 

Interestingly, for A > the ground-state energy is not 
a monotonic function of 77. It reaches its global minimum 
at some rj m i n which is smaller than ?7i_>.ii, i.e. the particles 
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Table 1: The values of the dimensionless layer separation 77 at 
phase transitions between structures I >• II >• III — > IV >■ V, 
for a given dielectric jump A. 



A 






f?IIT-»rV 


??IV^V 


0.2 





0.273 


0.619 


0.732 


0.4 





0.267 


0.607 


0.723 


0.6 





0.259 


0.596 


0.712 


0.8 





0.252 


0.587 


0.703 


1.0 





0.245 


0.578 


0.693 



form the structure-I Wigner monolayer at this point, see 
fig. |U Such a behavior differs substantially from the ho- 
mogeneous case A = 0, for which the energy grows mono- 
tonically. 
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Fig. 4: The dimensionless ground-state energy per particle, 
Ee/(Ne 2 y/2a), versus the dimensionless distance between the 
walls 77, for four values of the dielectric jump A = 0, 0.1, 0.4 
and 1. The symbols represent the transition points between 
structures: open circle for I —t II, open square for II — > III, 
filled diamond for III — > IV and filled triangle for IV — > V. 
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0.2 0.4 0.6 0.8 

Fig. 5: The aspect ratio 8 of structures I-III versus the di- 
mensionless distance between the walls r\. The dielectric jump 
A ranges between to 1, from left to right with step 0.1. 
5 = a/3 = 1.73205 ... for structure I, 1 < 5 < V3 within 
structure II and 8=1 for structure III. 




0.5 0.6 0.7 0.8 0.9 1 
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Fig. 6: Dependence of the angle ip for structure IV on the di- 
mensionless distance between the walls r\. The dielectric jump 
A ranges between to 1, from left to right with step 0.1. Each 
line ends at the transition point t?iv->v- 



At point 77i-j.11, the particle system starts to bifurcate to 
a Wigner bilayer (fj > 0) and at the same time structure 
I evolves into structure II with the varying aspect ratio 
1 < 5 < V3 (see fig. [SJ. Structure II ends up at the 
second-order transition point 77n_j.n1, where structure III 
with the fixed aspect ratio 5 = 1 determines the ground 
state energy. 

The second-order phase transition from structure III 
(ip = n/2) to structure IV, with a varying angle tp such 
that tan(y>/2) < 1, takes place at point 77m_j.1v. The plots 
of the order parameter tan(<p/2) versus the dimensionless 
inter wall distance are pictured for the dielectric jumps 
A = 0,0.1,..., 1 in fig. ED 

At the first-order transition point 77iv_>v, the angle ip 
skips from the value ~ 1.21 within structure IV, which is 



practically independent of A, to 7r/3 ~ 1.05 for structure 
V. An interesting question is whether also the layer sep- 
aration 77 exhibits a skip at 77iv_>v- The values of 771V-5.V 
for various A's in the last column of table Q] turn out to be 
the same for both structures IV and V at the transition 
point, within our number specification to 3 decimal digits. 

In the last fig. we plot the dimensionless distance 
of the Wigner sub-layer from the nearest wall, £, versus 
77, for four values of the dielectric jump A e [0, 1]. The 
dielectric homogeneous case A = with the trivial £ = 
0, i.e. particle collapse onto the wall surface, is shown 
for comparison. For the corresponding transition points 
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Fig. 7: Dimensionless distance of the Wigner sub-layer from 
the nearest wall, £, versus the dimensionless distance between 
the walls rj, for four values of the dielectric jump A = 0, 0.1, 0.4 
and 1. The dotted line has slope 1/2. For A = 0, we have the 
trivial result £ = 0, i.e. particles collapse onto the wall surface, 
as required by Earnshaw theorem. The horizontal dashed lines 
correspond to the asymptotic r\ — > oo values for the one-wall 
problem £*, given by (|13p . 



between the structures, see fig. |U If A > 0, the Wigner 
monolayer with structure I stays at z = up to 77 i_j.n and 
so £ = 77/2 for < 77 < Interestingly, the value of 

£ attained at f71_j.11 is always larger than the asymptotic 
77 — > 00 value £* for the one-wall problem, given by Eq. 
(fT3|) and depicted by horizontal dashed lines in fig. [Jj 
This is why, after bifurcation of the Wigner monolayer to 
the series of bilaycr structures II- V by increasing 77, the 
distance £ goes to £* from above. For r\ ~ 1.3, £ coincides 
with the asymptotic value £* to within 4 decimal digits. 

Conclusion. — We have extended the standard study 
of the homogeneous dielectric version of the ground state 
problem for particles with Coulomb pair interactions, con- 
strained between two equivalently charged wall, to a more 
realistic inhomogencous dielectric situation with repulsive 
electrostatic image charges. The analysis of the addi- 
tional Coulomb lattice sums due to the images of particles 
and surface charges was performed. All lattice sums are 
treated within the recent method [TSHU] which uses a bi- 
layer generalization of Misra functions ([TTj) to construct 
very quickly convergent series expansions for the ground 
state energy. The numerical results obtained show that 
inclusion of repulsive images has two fundamental effects. 
Firstly, the particles do not collapse onto the charged sur- 
faces of the walls, but create a Wigner monolayer or bilayer 
symmetrically with a nonzero spacing from the walls. The 
Wigner monolayer with the hexagonal structure I remains 
stable up to a certain (dimensionless) distance between 
the walls 771-j.ii, where it starts to bifurcate into the series 
of bilayer structures II-V. Secondly, the ground state en- 
ergy exhibits as the function of the distance between the 



walls a minimum at some 77 m ; n lying in the region of the 
structure-I Wigner monolayer. 

Although thermal fluctuations are supposed to destroy 
a classical Wigner crystal at very low temperatures, their 
fingerprints persist within a recent strong-coupling theory 
|24j to higher (room) temperatures, in agreement with MC 
simulations. Based on the present work, we plan to extend 
that strong-coupling study by including wall images. 
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